Libraries
Load data
I am using the AHA dataset for writing my code because it is the
simplest dataset I can work with. Simplest = significant number of
baseline samples without issue of repeated measures, balanced
distribution of classes, likely difference to be detected present.
p_ps = readRDS("/Users/aa370/Library/CloudStorage/Box-Box/project_davidlab/LAD_LAB_Personnel/Ammara_A/Projects/POMMS/R24_POMMS/data/processed/20221213_pomms_cleanmetadata_plant.rds")
p_ps
phyloseq-class experiment-level object
otu_table() OTU Table: [ 252 taxa and 384 samples ]
sample_data() Sample Data: [ 384 samples by 28 sample variables ]
tax_table() Taxonomy Table: [ 252 taxa by 8 taxonomic ranks ]
sample_data(p_ps)
Sample Data: [384 samples by 28 sample variables]:
c_ps = readRDS("/Users/aa370/Library/CloudStorage/Box-Box/project_davidlab/LAD_LAB_Personnel/Ammara_A/Projects/CHOICE/CHOICE_Trnl/CHOICE_20220912/ps_objects/choice_complete_20221205_clean.rds")
c_ps
phyloseq-class experiment-level object
otu_table() OTU Table: [ 571 taxa and 171 samples ]
sample_data() Sample Data: [ 171 samples by 30 sample variables ]
tax_table() Taxonomy Table: [ 571 taxa by 10 taxonomic ranks ]
Preprocess
Remove problematic samples
Taking out samples with less than or equal to 0 reads and keeping
only entry timepoint reduces samples from 384 down to 240.
p_ps = p_ps %>%
subset_samples( reads > 0) %>%
subset_samples(timepoint == "Entry")
p_ps
phyloseq-class experiment-level object
otu_table() OTU Table: [ 252 taxa and 240 samples ]
sample_data() Sample Data: [ 240 samples by 28 sample variables ]
tax_table() Taxonomy Table: [ 252 taxa by 8 taxonomic ranks ]
Class distribution in this dataset: Case = 193 Control = 47
sample_data(p_ps) %>%
as.data.frame() %>%
as_tibble() %>%
group_by(group) %>%
summarise(count = n())
Merge phyloseq
combined = merge_phyloseq(c_ps, p_ps)
combined
phyloseq-class experiment-level object
otu_table() OTU Table: [ 722 taxa and 411 samples ]
sample_data() Sample Data: [ 411 samples by 38 sample variables ]
tax_table() Taxonomy Table: [ 722 taxa by 10 taxonomic ranks ]
Aligning seq table
seqtab = combined@otu_table
seqtab_merged = dada2::collapseNoMismatch(seqtab)
combined = phyloseq(otu_table(seqtab_merged, taxa_are_rows=FALSE),
sample_data(combined@sam_data),
tax_table(combined@tax_table))
combined
phyloseq-class experiment-level object
otu_table() OTU Table: [ 684 taxa and 411 samples ]
sample_data() Sample Data: [ 411 samples by 38 sample variables ]
tax_table() Taxonomy Table: [ 684 taxa by 10 taxonomic ranks ]
#sample_data(combined)
Feature table
Zero variance features
No features have zero variance
nzv = features %>%
select(!sample) %>%
nearZeroVar(saveMetrics = TRUE)
nzv %>%
filter(nzv == TRUE)
rm(nzv)
Modeling
Loop
label = "group"
importance_df = data.frame()
cummulative_external = data.frame()
auc_list = c()
current_data = input %>%
select(-label_list[label_list != label]) %>%
select(-sample)
for(i in 1:20) {
# Data splitting
train_ind = caret::createDataPartition((current_data %>%
pull(label)), #split based on label column
p = 0.8,
list = FALSE,
times = 1)
train = current_data[train_ind,]
test = current_data[-train_ind,]
# Upsampling
train[[label]] = as.factor(train[[label]])
train = train %>%
as.data.frame()
formula_smote = label %>%
paste0("~.") %>%
as.formula()
train = DMwR::SMOTE(formula_smote, data = train)
# Training model
formula = label %>%
paste0("~.") %>%
as.formula()
fit_control = caret::trainControl(method = "CV", # 5-fold cross validation
number = 5,
summaryFunction = prSummary,
classProbs = T)
model = caret::train(formula, #using formula specified
data = train,
method = "rf",
trControl = fit_control,
verbose = FALSE,
metric = "AUC")
# Model predictions
predictions = predict(model, newdata = test, type= "prob")
predictions$obs = test %>%
pull(label) %>%
as.factor()
# Evaluate model
eval = MLeval::evalm(predictions, showplots = FALSE, silent = TRUE)
auc = eval$stdres$`Group1`[13,1]
auc_list = append(auc_list, auc)
#Extract feature importance
feature_importance = caret::varImp(model$finalModel) %>% # Ben calls this importances instead of final model but the rankings are exactly the same, only number value scale is different.
arrange(desc(Overall)) %>%
as.matrix() %>%
t() %>%
as.data.frame()
feature_importance$auc = auc
feature_importance$iteration = i
feature_importance$label = label
importance_df = rbind(importance_df, feature_importance)
#Predicting on external dataset
external_predictions = predict(model, newdata = external_test, type= "prob")
external_predictions$sample = external_test$sample
external_predictions$auc = auc
external_predictions$pred = predict(model, newdata = external_test)%>%
as.factor ()
external_predictions$iteration = i
cummulative_external = rbind(cummulative_external, external_predictions)
}
auc_list
[1] 0.61 0.61 0.48 0.76 0.74 0.32 0.82 0.46 0.68 0.70 0.78 0.68 0.70
[14] 0.70 0.56 0.57 0.42 0.68 0.65 0.57
importance_df
cummulative_external
NA
Plotting
AUC plot and t-test
data.frame(auc = auc_list)%>%
ggplot()+
geom_boxplot(aes(y = auc))+
theme_classic()

Save external predictions
predictions
write_csv(save, "/Users/aa370/Library/CloudStorage/Box-Box/project_davidlab/LAD_LAB_Personnel/Ammara_A/Projects/Machine_learning/LAD_lab_ML/code/pomms_to_choice_predictions.csv")
variable importances
write_csv(importance_df, "/Users/aa370/Library/CloudStorage/Box-Box/project_davidlab/LAD_LAB_Personnel/Ammara_A/Projects/Machine_learning/LAD_lab_ML/code/pomms_to_choice_feature_importance.csv")
LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKIyBMaWJyYXJpZXMKCmBgYHtyLCBpbmNsdWRlID0gRkFMU0V9CmxpYnJhcnkodGlkeXZlcnNlKQpsaWJyYXJ5KG1pY3JvYmlvbWUpCmxpYnJhcnkocGh5bG9zZXEpCmxpYnJhcnkoY2FyZXQpCmxpYnJhcnkocmFuZG9tRm9yZXN0KQpsaWJyYXJ5KE1CdXRpbHMpCmBgYAoKIyBMb2FkIGRhdGEKCkkgYW0gdXNpbmcgdGhlIEFIQSBkYXRhc2V0IGZvciB3cml0aW5nIG15IGNvZGUgYmVjYXVzZSBpdCBpcyB0aGUgc2ltcGxlc3QgZGF0YXNldCBJIGNhbiB3b3JrIHdpdGguClNpbXBsZXN0ID0gc2lnbmlmaWNhbnQgbnVtYmVyIG9mIGJhc2VsaW5lIHNhbXBsZXMgd2l0aG91dCBpc3N1ZSBvZiByZXBlYXRlZCBtZWFzdXJlcywgYmFsYW5jZWQgZGlzdHJpYnV0aW9uIG9mIGNsYXNzZXMsIGxpa2VseSBkaWZmZXJlbmNlIHRvIGJlIGRldGVjdGVkIHByZXNlbnQuCgoKYGBge3J9CnBfcHMgPSByZWFkUkRTKCIvVXNlcnMvYWEzNzAvTGlicmFyeS9DbG91ZFN0b3JhZ2UvQm94LUJveC9wcm9qZWN0X2RhdmlkbGFiL0xBRF9MQUJfUGVyc29ubmVsL0FtbWFyYV9BL1Byb2plY3RzL1BPTU1TL1IyNF9QT01NUy9kYXRhL3Byb2Nlc3NlZC8yMDIyMTIxM19wb21tc19jbGVhbm1ldGFkYXRhX3BsYW50LnJkcyIpCgpwX3BzCgpzYW1wbGVfZGF0YShwX3BzKQpgYGAKCmBgYHtyfQpjX3BzID0gcmVhZFJEUygiL1VzZXJzL2FhMzcwL0xpYnJhcnkvQ2xvdWRTdG9yYWdlL0JveC1Cb3gvcHJvamVjdF9kYXZpZGxhYi9MQURfTEFCX1BlcnNvbm5lbC9BbW1hcmFfQS9Qcm9qZWN0cy9DSE9JQ0UvQ0hPSUNFX1RybmwvQ0hPSUNFXzIwMjIwOTEyL3BzX29iamVjdHMvY2hvaWNlX2NvbXBsZXRlXzIwMjIxMjA1X2NsZWFuLnJkcyIpCgpjX3BzCmBgYAoKIyBQcmVwcm9jZXNzCgojIyBSZW1vdmUgcHJvYmxlbWF0aWMgc2FtcGxlcwoKVGFraW5nIG91dCBzYW1wbGVzIHdpdGggbGVzcyB0aGFuIG9yIGVxdWFsIHRvIDAgcmVhZHMgYW5kIGtlZXBpbmcgb25seSBlbnRyeSB0aW1lcG9pbnQgcmVkdWNlcyBzYW1wbGVzIGZyb20gMzg0IGRvd24gdG8gMjQwLgoKYGBge3J9CnBfcHMgPSBwX3BzICU+JQogIHN1YnNldF9zYW1wbGVzKCByZWFkcyA+IDApICU+JQogIHN1YnNldF9zYW1wbGVzKHRpbWVwb2ludCA9PSAiRW50cnkiKQoKcF9wcwpgYGAKCkNsYXNzIGRpc3RyaWJ1dGlvbiBpbiB0aGlzIGRhdGFzZXQ6CkNhc2UgPSAxOTMKQ29udHJvbCA9IDQ3CgpgYGB7cn0Kc2FtcGxlX2RhdGEocF9wcykgJT4lCiAgYXMuZGF0YS5mcmFtZSgpICU+JQogIGFzX3RpYmJsZSgpICU+JQogIGdyb3VwX2J5KGdyb3VwKSAlPiUKICBzdW1tYXJpc2UoY291bnQgPSBuKCkpCmBgYAojIyBNZXJnZSBwaHlsb3NlcQoKYGBge3J9CmNvbWJpbmVkID0gbWVyZ2VfcGh5bG9zZXEoY19wcywgcF9wcykKCmNvbWJpbmVkCmBgYAoKIyMjIEFsaWduaW5nIHNlcSB0YWJsZQoKYGBge3J9CnNlcXRhYiA9IGNvbWJpbmVkQG90dV90YWJsZQoKc2VxdGFiX21lcmdlZCA9IGRhZGEyOjpjb2xsYXBzZU5vTWlzbWF0Y2goc2VxdGFiKQoKY29tYmluZWQgPSBwaHlsb3NlcShvdHVfdGFibGUoc2VxdGFiX21lcmdlZCwgdGF4YV9hcmVfcm93cz1GQUxTRSksCiAgICAgICAgICAgICAgIHNhbXBsZV9kYXRhKGNvbWJpbmVkQHNhbV9kYXRhKSwgCiAgICAgICAgICAgICAgIHRheF90YWJsZShjb21iaW5lZEB0YXhfdGFibGUpKQoKY29tYmluZWQKCiNzYW1wbGVfZGF0YShjb21iaW5lZCkKYGBgCgoKIyMgRGF0YSB0cmFuc2Zvcm0KCmBgYHtyfQpvdHVfY2xyID0gYWJ1bmRhbmNlcyhjb21iaW5lZCwgImNsciIpICU+JSAKICB0KCkKCmNvbWJpbmVkID0gcGh5bG9zZXEob3R1X3RhYmxlKG90dV9jbHIsIHRheGFfYXJlX3Jvd3MgPSBGQUxTRSksCiAgICAgICAgICAgICAgICAgICBzYW1wbGVfZGF0YShjb21iaW5lZEBzYW1fZGF0YSksCiAgICAgICAgICAgICAgICAgICB0YXhfdGFibGUoY29tYmluZWRAdGF4X3RhYmxlKSkKCnJtKG90dV9jbHIpCmBgYAoKIyMgRmVhdHVyZSB0YWJsZQoKIyMjIEV4dHJhY3QgZGF0YQoKT3JkZXIgb2YgYXN2IGNvbG5hbWVzIGluIE9UVSB0YWJsZSBhbmQgdGF4IHRhYmxlIGlzIHRoZSBzYW1lLgoKYGBge3J9CmZlYXR1cmVzID0gY29tYmluZWRAb3R1X3RhYmxlICU+JQogIGFzLmRhdGEuZnJhbWUoKSAlPiUKICBhc190aWJibGUocm93bmFtZXMgPSBOQSkgJT4lCiAgcm93bmFtZXNfdG9fY29sdW1uKCJzYW1wbGUiKQoKZmVhdHVyZV9sYWJlbHMgPSBjKCJzYW1wbGUiKQoKZmVhdHVyZXMKYGBgCgojIyMgWmVybyB2YXJpYW5jZSBmZWF0dXJlcwoKTm8gZmVhdHVyZXMgaGF2ZSB6ZXJvIHZhcmlhbmNlCmBgYHtyfQpuenYgPSBmZWF0dXJlcyAlPiUKICBzZWxlY3QoIXNhbXBsZSkgJT4lCiAgbmVhclplcm9WYXIoc2F2ZU1ldHJpY3MgPSBUUlVFKQpgYGAKCmBgYHtyfQpuenYgJT4lCiAgZmlsdGVyKG56diA9PSBUUlVFKQpgYGAKCgpgYGB7cn0Kcm0obnp2KQpgYGAKCiMjIyBDb3JyZWxhdGVkIHByZWRpY3RvcnMKCmBgYHtyfQpmZWF0dXJlX2NvciA9IGZlYXR1cmVzICU+JQogc2VsZWN0KCFzYW1wbGUpICU+JQogY29yKCkgCgpzdW1tYXJ5KGZlYXR1cmVfY29yW3VwcGVyLnRyaShmZWF0dXJlX2NvcildKQpgYGAKCmBgYHtyfQpoaWdoX2NvciA9IGZpbmRDb3JyZWxhdGlvbihmZWF0dXJlX2NvciwgY3V0b2ZmID0gLjUwKSAjdXNpbmcgYSBtb3JlIHBlcm1pc3NpdmUgY3V0b2ZmIGJlY2F1c2UgYSBsb3Qgb2YgY29sdW1ucyBzdGlsbCBsZWZ0IGluIHRoZSBpbnB1dCB0YWJsZQoKZmlsdGVyZWRfZmVhdHVyZXMgPSBmZWF0dXJlc1ssLWhpZ2hfY29yXQpgYGAKCmBgYHtyfQpybShoaWdoX2NvcikKcm0oZmVhdHVyZV9jb3IpCmBgYAoKIyMjIExpbmVhciBkZXBlbmRlbmNpZXMKCmBgYHtyfQpmZWF0dXJlX2xkID0gZmlsdGVyZWRfZmVhdHVyZXMgJT4lCiAgc2VsZWN0KCFzYW1wbGUpICU+JQogIGZpbmRMaW5lYXJDb21ib3MoKQoKZmlsdGVyZWRfZmVhdHVyZXMgPSBmaWx0ZXJlZF9mZWF0dXJlc1ssIC0gZmVhdHVyZV9sZCRyZW1vdmVdCmBgYAoKYGBge3J9CnJtKGZlYXR1cmVfbGQpCmBgYAoKIyMgSW5wdXQgdGFibGUKYGBge3J9CnNhbXBsZV9sYWJlbHMgPSBjb21iaW5lZEBzYW1fZGF0YSU+JQogIGFzLmRhdGEuZnJhbWUoKSAlPiUKICBhc190aWJibGUocm93bmFtZXMgPSBOQSkgJT4lCiAgcm93bmFtZXNfdG9fY29sdW1uKCJzYW1wbGUiKSAlPiUKICBzZWxlY3Qoc2FtcGxlLCBncm91cCkgJT4lCiAgbGVmdF9qb2luKGZpbHRlcmVkX2ZlYXR1cmVzKQoKaW5wdXQgPSBzYW1wbGVfbGFiZWxzICU+JQogIGZpbHRlcihzdHJfZGV0ZWN0KHNhbXBsZSwgInYiKSkgJT4lICMgY2hhbmdlIGhlcmUKICBtdXRhdGUoZ3JvdXBfc2h1ZmZsZSA9IHNhbXBsZShncm91cCkpCgpleHRlcm5hbF90ZXN0ID0gc2FtcGxlX2xhYmVscyAlPiUKICBmaWx0ZXIoIXN0cl9kZXRlY3Qoc2FtcGxlLCAidiIpKSAlPiUgIyBjaGFuZ2UgaGVyZQogIHNlbGVjdCgtYyhncm91cCkpCmBgYAoKIyBNb2RlbGluZwoKIyMgQmFzZWxpbmUgcGVyZm9ybWFuY2UKCmBgYHtyfQpsYWJlbF9saXN0ID0gYygiZ3JvdXAiLCAiZ3JvdXBfc2h1ZmZsZSIpCnJlc3VsdHMgPSBsb29wX2xhYmVsKGl0ZXJhdGlvbnMgPTIwLCBpbnB1dCwgbGFiZWxfbGlzdCwgcGVyY2VudCA9IDAuOCwgdXBzYW1wbGluZyA9IFRSVUUpCgphdWNfZGYgPSByZXN1bHRzW1syXV0KCmltcG9ydGFuY2VfZGYgPSByZXN1bHRzW1sxXV0KYGBgCgpgYGB7cn0KYXVjX2RmICU+JQogIHBpdm90X2xvbmdlcihjb2xzID0gYygiZ3JvdXAiLCAiZ3JvdXBfc2h1ZmZsZSIpLCBuYW1lc190byA9ICJmZWF0dXJlIiwgdmFsdWVzX3RvID0gIkFVQyIpICU+JQogIGdncGxvdCgpKwogIGdlb21fYm94cGxvdChhZXMoeD1mZWF0dXJlLCB5ID0gQVVDKSkrCiAgdGhlbWVfY2xhc3NpYygpCgp0LnRlc3QoYXVjX2RmJGdyb3VwLCBhdWNfZGYkZ3JvdXBfc2h1ZmZsZSkKCmBgYAoKYGBge3J9CnJtKGF1Y19kZikKCnJtKGltcG9ydGFuY2VfZGYpCmBgYAoKIyMgTG9vcAoKYGBge3J9CmxhYmVsID0gImdyb3VwIgoKaW1wb3J0YW5jZV9kZiA9IGRhdGEuZnJhbWUoKQoKY3VtbXVsYXRpdmVfZXh0ZXJuYWwgPSBkYXRhLmZyYW1lKCkKCmF1Y19saXN0ID0gYygpCgpjdXJyZW50X2RhdGEgPSBpbnB1dCAlPiUKICAgICAgc2VsZWN0KC1sYWJlbF9saXN0W2xhYmVsX2xpc3QgIT0gbGFiZWxdKSAlPiUKICAgICAgc2VsZWN0KC1zYW1wbGUpCgpmb3IoaSBpbiAxOjIwKSB7CiAgIyBEYXRhIHNwbGl0dGluZwogICAgdHJhaW5faW5kID0gY2FyZXQ6OmNyZWF0ZURhdGFQYXJ0aXRpb24oKGN1cnJlbnRfZGF0YSAlPiUKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHB1bGwobGFiZWwpKSwgI3NwbGl0IGJhc2VkIG9uIGxhYmVsIGNvbHVtbgogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgcCA9IDAuOCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGxpc3QgPSBGQUxTRSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpbWVzID0gMSkKICAgIHRyYWluID0gY3VycmVudF9kYXRhW3RyYWluX2luZCxdCiAgICB0ZXN0ID0gY3VycmVudF9kYXRhWy10cmFpbl9pbmQsXQogICAgCiAgICAjIFVwc2FtcGxpbmcKICAgIHRyYWluW1tsYWJlbF1dID0gYXMuZmFjdG9yKHRyYWluW1tsYWJlbF1dKQogICAgdHJhaW4gPSB0cmFpbiAlPiUKICAgICAgICBhcy5kYXRhLmZyYW1lKCkKICAgICAgCiAgICBmb3JtdWxhX3Ntb3RlID0gbGFiZWwgJT4lCiAgICAgIHBhc3RlMCgifi4iKSAlPiUKICAgICAgYXMuZm9ybXVsYSgpCiAgICAgIAogICAgdHJhaW4gPSBETXdSOjpTTU9URShmb3JtdWxhX3Ntb3RlLCBkYXRhICA9IHRyYWluKQogICAgCiAgICAjIFRyYWluaW5nIG1vZGVsCiAgICAKICAgIGZvcm11bGEgPSBsYWJlbCAlPiUKICAgICAgcGFzdGUwKCJ+LiIpICU+JQogICAgICBhcy5mb3JtdWxhKCkKICAgIAogICAgZml0X2NvbnRyb2wgPSBjYXJldDo6dHJhaW5Db250cm9sKG1ldGhvZCA9ICJDViIsICMgNS1mb2xkIGNyb3NzIHZhbGlkYXRpb24KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbnVtYmVyID0gNSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc3VtbWFyeUZ1bmN0aW9uID0gcHJTdW1tYXJ5LAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjbGFzc1Byb2JzID0gVCkKICAgIAogICAgbW9kZWwgPSBjYXJldDo6dHJhaW4oZm9ybXVsYSwgI3VzaW5nIGZvcm11bGEgc3BlY2lmaWVkIAogICAgICAgICAgICAgIGRhdGEgPSB0cmFpbiwKICAgICAgICAgICAgICAgICAgbWV0aG9kID0gInJmIiwKICAgICAgICAgICAgICAgICAgdHJDb250cm9sID0gZml0X2NvbnRyb2wsCiAgICAgICAgICAgICAgICAgIHZlcmJvc2UgPSBGQUxTRSwgCiAgICAgICAgICAgICAgICAgIG1ldHJpYyA9ICJBVUMiKQogICAgCiAgICAjIE1vZGVsIHByZWRpY3Rpb25zCiAgICAKICAgIHByZWRpY3Rpb25zID0gcHJlZGljdChtb2RlbCwgbmV3ZGF0YSA9IHRlc3QsIHR5cGU9ICJwcm9iIikKICAgIAogICAgcHJlZGljdGlvbnMkb2JzID0gdGVzdCAlPiUKICAgICAgcHVsbChsYWJlbCkgJT4lCiAgICAgIGFzLmZhY3RvcigpIAogICAgCiAgICAjIEV2YWx1YXRlIG1vZGVsCiAgICAKICAgIGV2YWwgPSBNTGV2YWw6OmV2YWxtKHByZWRpY3Rpb25zLCBzaG93cGxvdHMgPSBGQUxTRSwgc2lsZW50ID0gVFJVRSkKICAgIAogICAgYXVjID0gZXZhbCRzdGRyZXMkYEdyb3VwMWBbMTMsMV0KICAgIAogICAgYXVjX2xpc3QgPSBhcHBlbmQoYXVjX2xpc3QsIGF1YykKICAgIAogICAgI0V4dHJhY3QgZmVhdHVyZSBpbXBvcnRhbmNlCgogICAgZmVhdHVyZV9pbXBvcnRhbmNlID0gY2FyZXQ6OnZhckltcChtb2RlbCRmaW5hbE1vZGVsKSAlPiUgIyBCZW4gY2FsbHMgdGhpcyBpbXBvcnRhbmNlcyBpbnN0ZWFkIG9mIGZpbmFsIG1vZGVsIGJ1dCB0aGUgcmFua2luZ3MgYXJlIGV4YWN0bHkgdGhlIHNhbWUsIG9ubHkgbnVtYmVyIHZhbHVlIHNjYWxlIGlzIGRpZmZlcmVudC4KICAgICAgYXJyYW5nZShkZXNjKE92ZXJhbGwpKSAlPiUKICAgICAgYXMubWF0cml4KCkgJT4lCiAgICAgIHQoKSAlPiUKICAgICAgYXMuZGF0YS5mcmFtZSgpCiAgICAKICAgIGZlYXR1cmVfaW1wb3J0YW5jZSRhdWMgPSBhdWMKICAgIAogICAgZmVhdHVyZV9pbXBvcnRhbmNlJGl0ZXJhdGlvbiA9IGkKICAgIAogICAgZmVhdHVyZV9pbXBvcnRhbmNlJGxhYmVsID0gbGFiZWwKICAgIAogICAgaW1wb3J0YW5jZV9kZiA9IHJiaW5kKGltcG9ydGFuY2VfZGYsIGZlYXR1cmVfaW1wb3J0YW5jZSkKICAgIAogICAgI1ByZWRpY3Rpbmcgb24gZXh0ZXJuYWwgZGF0YXNldAogICAgCiAgICBleHRlcm5hbF9wcmVkaWN0aW9ucyA9ICBwcmVkaWN0KG1vZGVsLCBuZXdkYXRhID0gZXh0ZXJuYWxfdGVzdCwgdHlwZT0gInByb2IiKQogICAgCiAgICBleHRlcm5hbF9wcmVkaWN0aW9ucyRzYW1wbGUgPSBleHRlcm5hbF90ZXN0JHNhbXBsZQogICAgCiAgICBleHRlcm5hbF9wcmVkaWN0aW9ucyRhdWMgPSBhdWMKICAgIAogICAgZXh0ZXJuYWxfcHJlZGljdGlvbnMkcHJlZCA9IHByZWRpY3QobW9kZWwsIG5ld2RhdGEgPSBleHRlcm5hbF90ZXN0KSU+JQogICAgICBhcy5mYWN0b3IgKCkKICAgIAogICAgZXh0ZXJuYWxfcHJlZGljdGlvbnMkaXRlcmF0aW9uID0gaQogICAgCiAgICBjdW1tdWxhdGl2ZV9leHRlcm5hbCA9IHJiaW5kKGN1bW11bGF0aXZlX2V4dGVybmFsLCBleHRlcm5hbF9wcmVkaWN0aW9ucykKICAgIAp9CgphdWNfbGlzdAoKaW1wb3J0YW5jZV9kZgoKY3VtbXVsYXRpdmVfZXh0ZXJuYWwgIAogIApgYGAKCiMgUGxvdHRpbmcKCiMjIEFVQyBwbG90IGFuZCB0LXRlc3QKYGBge3J9CmRhdGEuZnJhbWUoYXVjID0gYXVjX2xpc3QpJT4lCiAgZ2dwbG90KCkrCiAgZ2VvbV9ib3hwbG90KGFlcyh5ID0gYXVjKSkrCiAgdGhlbWVfY2xhc3NpYygpCmBgYAoKIyMgU2F2ZSBleHRlcm5hbCBwcmVkaWN0aW9ucwoKIyMjIHByZWRpY3Rpb25zCmBgYHtyfQpzYW1kZiA9IGNfcHNAc2FtX2RhdGElPiUKICBhcy5kYXRhLmZyYW1lKCkgJT4lCiAgYXNfdGliYmxlKHJvd25hbWVzID0gTkEpICU+JQogIHJvd25hbWVzX3RvX2NvbHVtbigic2FtcGxlIikgJT4lCiAgc2VsZWN0KHNhbXBsZSwgaWQsIHRyZWF0bWVudCwgdHJ1ZV93ZWVrLCB0aW1lcG9pbnQpIAoKc2F2ZSA9IGN1bW11bGF0aXZlX2V4dGVybmFsICU+JQogICNmaWx0ZXIoYXVjID09IG1heCkgJT4lCiAgbGVmdF9qb2luKHNhbWRmKSAlPiUKICBhcnJhbmdlKGl0ZXJhdGlvbiwgdHJlYXRtZW50LCBpZCwgdHJ1ZV93ZWVrKQpgYGAKCmBgYHtyfQp3cml0ZV9jc3Yoc2F2ZSwgIi9Vc2Vycy9hYTM3MC9MaWJyYXJ5L0Nsb3VkU3RvcmFnZS9Cb3gtQm94L3Byb2plY3RfZGF2aWRsYWIvTEFEX0xBQl9QZXJzb25uZWwvQW1tYXJhX0EvUHJvamVjdHMvTWFjaGluZV9sZWFybmluZy9MQURfbGFiX01ML2NvZGUvcG9tbXNfdG9fY2hvaWNlX3ByZWRpY3Rpb25zLmNzdiIpCmBgYAoKCiMjIyB2YXJpYWJsZSBpbXBvcnRhbmNlcwoKYGBge3J9CndyaXRlX2NzdihpbXBvcnRhbmNlX2RmLCAiL1VzZXJzL2FhMzcwL0xpYnJhcnkvQ2xvdWRTdG9yYWdlL0JveC1Cb3gvcHJvamVjdF9kYXZpZGxhYi9MQURfTEFCX1BlcnNvbm5lbC9BbW1hcmFfQS9Qcm9qZWN0cy9NYWNoaW5lX2xlYXJuaW5nL0xBRF9sYWJfTUwvY29kZS9wb21tc190b19jaG9pY2VfZmVhdHVyZV9pbXBvcnRhbmNlLmNzdiIpCmBgYAo=